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. A version of the Gaussian self-consistent (GSC) method, which avoids the use of the Edwards' 

virial expansion, is presented. Instead, the mean energy is evaluated directly via a convolution of 
the attractive part of the pair-wise non-bonded potential with the Gaussian trial radial distribution 
function. The hard sphere repulsion is taken into account via a suitably generalised Carnahan- 
QQ ' Starling term. Comparison of the mean-squared inter-monomer distances and radius of gyration, 

as well as of the mean energy, between the results from the GSC calculations and Monte Carlo (MC) 
simulation in continuous space are made across the coil-to-globule transition for isolated ring, open 
I and star homopolymers of varied lengths and flexibility. Importantly, both techniques utilise the 

same polymer model so that the data points could be directly superimposed. A surprisingly good 
overall agreement is found between these GSC and MC results. Caveats of the Gaussian technique 
\ and ways for going beyond it are also discussed. 
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PACS numbers: 36.20.-r, 36.20.Ey, 61.25.Hq 



I. INTRODUCTION 

> 

' The Gaussian self-consistent (GSC) method represents a quite general, albeit relatively simple, theoretical frame- 
^sj , work for studying the equilibrium, dynamics and kinetics of conformational changes in polymer solutions. One of 
C ■ ■ its most attractive features is that it can be-|applied to virtually any type of heteropolymerEiB, with any specific 
interactijOn termsJnvolving e.g. chain stifFneasQcl or chargeaS. The connectivity of the iiain can be arbitrary: an opieii 
polymerQ, a ringj, a staiM, or a dendrimenlEl, whether in the limit of a single chaintiJ or at finite concentrationsllj. 
Moreover, not only the radius of gyration of the polymer or the total density can be computed, but much finer details 
of the conformational structure, such as the mean-squared distances between monomers or radial distributions at a 
given time, are available. Such a versatility naturally comes at the cost of certain limitations and inaccuracies intrinsic 
'Y' , to the technique. However, these are well known and manageable, although often largely overstated in the popular 
view. Some of these, such as the Gaussian shape of the radial distribution function (RDF), are qpte-, unavoidable at 
, the level of the GSC theory, but others, such as the traditional use of the Edwards' virial expansior£3iiJ, can indeed be 
O ■ surpassed as we shall demonstrate in the current work. This would permit us for the first time to compare observables 
between the predictions of the GSC technique and Monte Carlo (MC) simulation in continuous space based on the 
same model for given values of thermodynamic parameters throughout the range of the coil-to-globule transition. 

Historically, perhaps the first of the equilibrium versions of the jQaussian method in apjplication to a single ho- 
mopolymer chain Wre independently proposed by S. Edwards et atj and J. des Cloizeauxllj. The approach of the 
' work by EdwardsEj relied on the following idea. One would like to approximate the coordinates of the real ho- 
mopolymer chain, Xj, governed by the Gibbs distribution with the exact Hamiltonian H, via a 'trial' Gaussian chain 
with the monomer positions xf'\ governed by a simpler trial Hamiltonian H'^^^ By assuming that the quantities 

AXi — Xi — Xf-' and AH = H — are small and hence by applying the Taylor expansion to the equilibrium 
Gibbs averages, in the first order, one requires that the mean-squared radius of gyration obeys, 

{Rl)o = {Rf^\, or equivalently «) 2i/(°))o - 2)o (i/^"^ = 0. (1) 
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Such a theory predicts the correct mean-field swelHng exponent i/p = 3/5 for the repulsive coil, only missing the 
subtle renormalisation group correction. It also gives appropriately vg — 1/2 for the ideal coil, and Ug = 1/3 for the 
^Jo^jl^j^^oreover, this approach can be further extended to kinetics by considering time dependent analogues of Eq. 

A seemingly more refined theory proposed by des Cloizeauxlli uses — 1 (where N is the degree of polymerisation) 
variational variables xlj'\ which are the normal modes of the monomer positions for a ring chain xf \ This technique 
is based on the standard Gibbs-Bogoliubov variational principle with a quadratic diagonal trial Hamiltonian H^^\ 
The equilibrium corresponds to the minimum of the trial free energy. Atrial = A'^°^ + {H- H<-°^)o, where A^°^ is the 
free energy associated with H^^\ Expressed in the Edwards' form the extremum conditions for the trial free energy 
are precisely the following A'^ — 1 equations, 

(X2)o = (X(")2)o, or equivalently (X^ ^ffW)^ - (X^ (ijW)o = 0. (2) 

Unfortunately, this theory, although more accurate for polymers around the theta-point and indeed for the globule, 
is known to have a deficiency in that it predicts the exponent = 2/3 instead of the mean-field value of i>f = 3/5 for 
the swollen coil in the thermodynamic limit N oo. For this reason, some suspicion arose concerning the validity of 
the theory by des Cloizeaux, hampering further efforts in this direction. It should be emphasised, however, that for 
finite values of A^ the resulting theory does give quite reasonable numerical results fairly close to = 3/5. In fact, 
the numerically fitted value of the swelling exponent tends to be somewhere in between the two above theoretical 
predictions. In a way, this is not an unusual situation, even foij-the more elaborate integral equations theories of 
molecular fluids used in attempts to deduce the Flory exponentoEj. More importantly, this well known deficiency of 
the theory by des Cloizeaux becomes irrelevant when considering the theta-point or the poor solvent conditions. 

Therefore, activity o]xJuxt]it.'r /.^X|t.ci3ding and improving the GSC theory has persisted,. notably in Milan by G. Allegra 
and F.|-jGaMB«©li et a£il-Du'trE3nJ'E2l, in Saclay by H. Orland and T. jGarel et aOlla'El, as well as in our Group in 
DublinM^yiia, and some collaborations between these have emergedQ. Furthermore, extensions of the GSC theory 
to kinetics, as well as to essentially any polymeric system, have been also achieved. For the former end, one proceeds 
from the Langevin equation and approximates the exact stochastic ensemble via a linear trial one. We refer the reader 
to our papers in Refs. 0,|ll|,^ for the ultimate formulation of the GSC theory in terms of the real space variables. Such a 
version of the GSC method avoids the limitations of the earlier normal modes formulations since it can distinguish the 
frustrated states of hetcropolymers with spontaneously broken kinematic symmetries associated with the connectivity 
structure. 

Despite the relative success of the most recent form of the GSC method, the theoretical situation is not as yet fully 
satisfactory and further efforts along these lines are required. Firstly, one would like to be able to use the actual 
molecular interaction parameters rather than the coefficients in the Edwards' virial expansion, which are somewhat 
obscure in their meaning. Secondly, the convergence of such a virial expansion in the dense liquid globular state is 
rather problematic, being at best that of an asymptotic series. Thirdly, despite the conventional view, the three-body 
term is unable by itself to fully withstand a catastrophic pair-wise collapse of monomers onto each other in the 
attractive regime for hetcropolymers. In Refs. an ad hoc self-interaction energy term has been added to tackle 
this problem, but it is merely a tweak of the model. Finally, the effective three-body and higher order virial terms 
make the current theory numerically inefficient due to a scaling law involving a high power in A for the time expense 
per step, despite the fact that the original model only contains two-body interactions. 

Therefore, in the current work we shall rid ourselves of relying on the virial representation of the Hamiltonian 
entirely. Indeed, the exact molecular two-body potential can be averaged over the Gaussian RDF directly. However, 
the Gaussian function does not vanish at the origin, yielding a divergence from the repulsive part of the Lennard- Jones 
pair- wise potential. One simple alternative is to use the latter potential with a hard core part, as in Ref. in which 
case one has to cut off the integration at the hard sphere diameter, or more generally, within the excluded volume 
area. 

To include the effect of the hard sphere repulsion one can in principle use any of the standard reference techniques 
developed for the haxd-sphere liquids, e.g. the Percus-Yevick expression for RDF, or the Carnahan-Starling (CS) 
free energy formulaE^Ej. We prefer to use the latter due to its practical simplicity and accuracy when compared to 
molecular simulations. The Carnahan-Starling equation, which was originally derived for a single component hard 
sphere liquid, has been naturally extended to mixtures. It is quite popular in the liquid matter literature at present 
(see e.g. Ref. Its extension to polymeric fluids is somewhat less obvious though. Here we propose to express 
the partial packing coefficient for each of the monomers via an integral of the Gaussian trial RDF over the excluded 
volume area. The resulting GSC equations involve well tractable, albeit somewhat more complicated than before, 
expressions, since they now depend on the exact shapes of the molecular potentials. 

Although we shall present the resulting technique in its most general form, which is applicable to any type of poly- 
meric system, our numerical analysis in this work will be restricted by homopolymers of three different architectures 
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— ring, open and star — across the coil-to-globule transition. This is done to enable us to make a comparison of the 
results from the present GSC theory with the MC data obtained in our previous paper in Ref. 



II. MODEL 

The current coarse-grained polymer model is based on the following Hamiltonian (energy functional)iiflB@ 
terms of the monomer coordinates, X.^: 



^ Kij (X» - Xj)2 + — ^ Aijfe(Xi + Xfe - 2Xj 



\ J2 (f^S^'^dX. - X,|) + [/(;°"\|X, - X,|)) . (3) 
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Here the first term represents the connectivity structure of the polymer with harmonic springs of a given strength 
Kij, introduced between any pair of connected monomers (which is denoted by i ~ j). The second term represents 
the bending energy penalty given by the square of the local curvature with a characteristic stiffness Xijk between any 
three consecutively connected monomers (which is denoted by i « j « fc) in the form of the Kratky-Harris-Hearst. 
Below we shall prefer to rewrite the first two terms in the following equivalent form. 

Finally, the third and fourth terms represent pair-wise non-bonded interactions between monomers such as the van 
der Waals and Coulomb forces. We can adopt the Lennard-Jones form of the former potential, 

^ (0) , (0) 

hoo, r < rl ' + 



where there is also a hard core part with the monomer radii rl"' , and where U^j"^ are the dimensionless strengths of 
the interactions. The Coulomb interaction potential similarly is, 

C/r^M = 1. l/f eM~r/ln), h ^ ^^|f^, (6) 

where Id is the Debye screening length, eo is the dielectric permittivity of vacuum, Qi = Qq Qi arc the charges, and 
Ib is the Bjerrum length. 



III. METHOD 



A. Equations of the GSC method 



The main objects in the GSC method are the mean-squared distances between monomers, 

%(t)^i^(X,(t)-X,(Ory (7) 

Note that our convention includes the factor of 1 /3 here and hence in the later definition of the mean-squared radius 
of gyration according to the tradition. This allows us to rid of such factors from the radial distribution function 
(RDF) and numerous averages over it. 

The GSC method is based on replacing the stochastic ensemble for Xi with the exact Hamiltonian in the Langevin 
equation of motion onto the trial ensemble xf^\t) with a trial Hamiltonian H^^'^{t). The latter is taken as a generic 
quadratic form with the matrix coefficients, which are called the time-dependent effective potentials. 
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F(°)[x(0] = i^y.,(t)x,(t)x,(i). (8) 

ij 

Then one requires that the inter-monomer correlations satisfy the condition, 

{X,{t)X,{t)),^{xf\t)xf\t)),, (9) 

which means that the trial ensemble well approximates the exact one as far as monomer correlations are concerned. 
This procedure yields expressions for the effective potentials via the instantaneous mean energy, 

and hence the mean-squared distances themselves satisfy the self-consistent equationtS. These, in the absence of the 
hydrodynamic interaction, are simply, 



2 dt ' 3 ' " V dV,k{t) dV.kit) ) 



(11) 



Here C,h is the friction coefficient of a monomer, and the instantaneous free energy has the same functional expression 
via the instantaneous 2?ij(i) as it has at equilibrium. Extension to the preaveraged hydrodynamic approximation is 
quite straightforward also and it is discussed in Ref. ^ 

The stationary limit of these equations produces the equations for the minimum of the free energy, which are the 
same as those derived from the Gibbs-Bogoliubov variational principle with a generic quadratic trial Hamiltonian. 
Although in this paper we shall only be concerned with the equilibrium properties, the numerical solution of Eq. jll]), 
applied until the stationary limit is reached, presents by far the most efficient techniq|Ue for finding the global free 
energy minimum. This, based on the fifth order adaptive step Runge-Kutta integratom, was used for obtaining the 
results from the GSC technique in this work. 

It should be also noted that the systems studied here possess a large number of kinematic symmetries for TDij , and 
hence for Vij , matrices coming from their symmetricity and from the equivalence of any monomer in a ring, or any arm 
in a star homopolymers. Thus, the computational expenses per step in our calculations are of order tc ^ N F, where 
N/2 ^ F < is the total number of independent elements in the matrix 2?^. These symmetries significantly 

reduce the computational times compared to MC for an equivalent system, where such symmetries only appear in the 
observables after averaging over the statistical ensemble. For comparison, the computational expenses per step in MC 
are of order tc ^ N At S, where At ~ iV^ is the number of MC steps needed to ensure a good statistical independence 
between measurements, and where S is the number of measurements needed for sampling of observables. Typical 
values of S should be of order of 10"* — 10^ for a good accuracy in the present casecfl. Moreover, kinetic iteration of 
the GSC equations towards the equilibrium is also significantly faster than the equivalent equilibration procedure in 
the MC. For example, foif--a ring polymer of = 150 units in the good solvent, J/^"-' = 1, the GSC method turns out 
to be « 1000 times fasterB than MC for gaining the same data, whereas for an open chain of the same length, which 
possesses much fewer kinematic symmetries, GSC is « 200 times faster than MC. 



B. Hard sphere contribution 

Carnahan and Starling have devised a simple but rather accurate equation of state for hard sphere liquidB 
terms of the packing coefficient 77, yielding the following free energy. 



ksTN ' (l-r;)2 ' " ^ V ' ^ ^ 

where and V are the volumes of the hard sphere and of the whole system respectively. This is obtained from an 
interpolation formula for the virial coefficients based on several of them known exactly (see more detailed discussions 
in e.g. Ref. ||). 

To extend these ideas to polymeric fluids, we would have to distinguish the packing coefficients for individual 
monomers. In Eq. (|l2|) we can write N = the volume of the sphere is equal to the 1/8-th of its excluded 

volume; and the inverse volume is equal to the RDF of the ideal reference system, 1/V = gl^^^^i- Thus, we can similarly 
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express the packing coefRcient for the z-th monomer as a sum over all other monomers of the 1/8-th of the integral 
over the excluded volume of the monomer i of the RDF for the Gaussian reference system, 

yielding finally, 



\r\<T 



One may note that the Gaussian distribution does not possess a well defined volume beyond which it vanishes. Thus, 
to account for this we can, in principle, include a multiplicative parameter a, which should be once and for all chosen 
to match the data best, but this should, in any case, be fairly close to the unity. Therefore, the total hard sphere 
contribution will be, 

^-NFi^'^\,) (15) 



C. Free energy in the GSC method 

The total mean energy includes both the bonded and the pair-wise non-bonded interactions via. 



2?.. + dv ^ (r) Pl ^^ (r) + U^r' (r)] . (16) 



Note that here we integrate only beyond the excluded volume as the hard sphere contribution will be included explicitly 
via the Carnahan-Starling term. 

Given that the conformational Gaussian entropy has been calculated by us in Appendix B of Rcf. |l|, the total free 
energy can be summarised as follows, 

A = £-TS, £ = =£(f'°"'i) 5 = 5(3'^"), (17) 

where its various terms are given by, 

c{bond) 

=3EC/J°"''^I?.„ (18) 



i<j 



j:uS^E(^'^bJ^,], (19) 



ksT 

kgT " ^ JO) (0) 
i<j ' i ^ ' j 

^ ^lndeti?(^-i), (21) 



E4f^^^^™'[2/^.'(^^°^+-r)//z.], (20) 



kB 2 



and the arguments of these functions are defined as. 



Rij - 



.(0) , „(o) ' 



3 ij 



(22) 



r,, = aJ2F^''Hy^3h (23) 

(24) 
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Finally, the functions of yij are, 



dx exp ( — 



2y2 



8y 



1 2 



y \ Jl 



1 1 



/27^ exp (-2!^) (y - + 3y5 + 300y^ - 210y^) ~ erfc (^) (l + SlSy^ 



1 /2 



exp 



X da; exp 



^ exp(— fcx) 



1 



2y2 



1/2 2fc 



y \ TT TT 



■ exp 



(1 + ky 
2y2 



2\2 



erfc 



1 + ky"' 



y 



Here we have used the standard definitions for the error functions, 



(26) 



(27) 



(28) 



erf(z) = —= / dx exp(-a;^), erfc(z) = 1 - erf(z). (29) 
VTT Jo 

Note also that for a large z, numerically, one has to use the truncated asymptotic expansion in calculating E^'^°^^ to 
avoid divergences, 



erfc(z) exp(z2) 



1 



1 



1 3 

2^2 + 



15 

8^6 



(30) 



IV. RESULTS 

Here we shall restrict ourselves by the case of homopolymers, so that in Eq. (|^) all non-zero bonded interaction 
constants are equal: = k, Ay^ = A, as well as all non-bonded interaction parameters are identical in Eq. (j^): 
uf^'^ — U^^\ rf^"^ = d/2. We also choose the hard sphere diameter d equal to the length £ defined in Eq. (||), as in 

Ref. Moreover, henceforth we shall use the mean energy expressed in units of fcsT and the mean-squared 

distances Vij and the mean-squared radius of gyration 3TZ^ expressed in units of i'^. 

Firstly, to understand the infiuence of the adopted Carnahan-Starling term, we shall look at the case of a ring 
homopolymer with varied values of the spring constant in the good athermal solvent, U^'^^ = 0. In Tab. | we compare 
values of the mean energy between MC and GSC with two choices of the multiplicative parameter a in Eq. ( [2^ ) 

for the packing coefficient, namely a — 0.9 and a — 1.0. Likewise, in Tab. ^we present the data for the mean-squared 
radius of gyration 37?,g. One can see that the results from the GSC theory simply coincide with those from MC as 
K — > 0. Yet, the agreement in the energy is somewhat better for the theory with a = 0.9 than for that with the naive 
choice a = 1.0, and that is how exactly the particular value a = 0.9 has been chosen by us. A typical relative error in 
the energy is less than 1.5 percent for a — 0.9 and is under 4 percent for a = 1.0. The agreement for the mean-squared 
radius of gyration is somewhat less impressive — the relative error steadily increases with n, being better for a — 0.9, 
but never exceeding a few dozen percent. 

One can comment on the reason why a value a < I produces a somewhat better agreement with MC. To facilitate 
this discussion, RDFs from the GSC and MC techniques are exhibited in Fig. ^ In Ref. |2^ we have discussed the role 
of the 'correlation hole' at small separations beyond the excluded volume area, which is also evident for r ^ 15 in the 
main part of Fig^U As this feature is absent in the Gaussian theory (see Eq. (p^)), the resulting packing coefficient 
values in Eq. (p2) are overestimated. Thus, reducing the parameter a permits us to effectively lower the packing 
coefficient and this can be done only once because the correlation hole effect is expressed in terms of the function of 
the dimensionless variables, 

(31) 
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Next, we would like to analyse the dependence on the degree of polymerisation, N, in the good solvent, U^^^ = 0, 
(the second block in Tabs. |, The mean energy in this case comes entirely from the bonded interactions and 
hence the relative energy deviation of GSC from MC does not increase with N. This means simply that the values 
of match with those from MC quite well. As for the radius of gyration, the disagreement increases slowly with 

N, i.e. the mean-squared distances I?i,i+fe are increasingly more overestimated by GSC as compared to MC for large 
k. Based on these results, we can determine the swelling exponent of the ring coil by fitting 37?.^ via Eq., 

37^^ = N^" (32) 

using our data in the range = 50 — 500. The resulting prefactor and the exponent ly in all three cases are, 

( 0.223 ±0.005 MC 
6^ = <^ 0.180 ±0.003 GSC, a = 0.9 (33) 



0.223 ±0.005 


MC 






0.180 ±0.003 


GSC, 


a — 


0.9 


0.184± 0.003 


GSC, 


a — 


1.0 


0.610 ±0.006 


MC 






0.648 ±0.002 


GSC, 


a = 


0.9 


0.651 ±0.001 


GSC, 


a — 


1.0 



(34) 

^ 0.651 ±0.001 GSC, a= 1.0 

Thus, the increase in the deviation of 3TZg in the GSC method from MC indeed comes entirely from the overestimation 
of the swelling exponent. This is believed to reach the value v = 2/3 — 0.666. . . asymptotically. However, in the 
present range of N both-QSC and MC give higher apparent exponent values than the most accurate renormalisation 
group result up-to-data53, v = 0.5882 ± 0.0011. The GSC results for v are only slightly sensitive on the value of a, 

(2) 

being in a relatively small overestimation over the MC results. This is related to the Gaussian shape of g^j in the 
GSC theory, whereas a stretched exponential tail, exp{—Bf^), (see the range r > 40 in Fig. ^ contributes most to 
3TZg in MC. We may note also that the GSC estimate for v is fairly close to the result f = 0.635 from the integral 
equations approach based on a complicated closure of the Born-Green- Yvon hierarchy in Ref. pO| . 

When considering open flexible homopolymers (k = 1, A = 0) in the good solvent (C/*^*^^ — 0), the general behaviour 
of and 37?.^ in the third block of Tabs. |, || is very similar to that of a ring: a good and almost A^-independent 

agreement for and a slow increase with A^ in the relative error for 37?.^. In Fig. || we plot the mean-squared 

distances from the end monomer I?ofc vs the chain index k. Up to fc ~ 10, the agreement of both GSC curves with 
the MC data is nearly perfect, whereas both GSC curves increasingly overestimate the MC data for larger k. This 
is consistent with the overestimation of 37?.^ by the GSC method, dominated by large k contributions. Notably, the 
effect of changing the parameter a is rather weak on this scale. 

Further, let us investigate the effect of increasing the chain stiffness A from a value corresponding to a fairly flexible 
ring A = 1 (see Fig. Band the fourth block in Tabs. |, to that of a semi-stiff ring A = 5 (see Fig. ^ and the 
fifth block in Tabs. ITy). The fairly flexible case gives the energies in a very good agreement with the MC data, 
albeit the theory with a ~ 0.9 is somewhat less accurate accidentally. However, the agreement of 37?.^ in the GSC 
theory with MC is even better here than for the corresponding flexible coil. The semi-flexible case also gives the 
energies in a good agreement with the MC data, while STZ^ tends to be underestimated by the GSC method. Plots 
of the mean-squared distances I?o/c vs the chain index k in Figs. ^, ^ match nearly perfectly up to A; ^ 10 as well, 
diverging for larger k. The MC curve is lower (higher) than the GSC curves for A = 1 (A = 5) in accord with the 
tables data. Thus, overall, the increase of 3TZg with the chain stiffness A is more rapid in MC than in the GSC theory. 
Indeed, conformations of a stiff chain in MC become those of a rigid ring (or rod for an open chain) with increasing A. 
However, the Carnahan-Starling equation was deduced in the assumption of a total 3-d isotropy, when its influence 
would be weaker than in an effective 1-d projection. 

Next, let us bring our attention to the effect of changing the topology of the chain. Thus, in the sixth block of 
Tabs. |, |l| we present and 3TZ^ for the flexible (A = 0) stars with the arm length N/ f = 50 in the good solvent 

(C/(o) ^ 0). The relative errors m and 37^2 increase with the number of arms / steadily, again the energy values 

being more close between MC and GSC. The mean-squared distances Dqi- from the core monomer for the largest star 
with f — 12 arms are plotted in Fig. ^. Clearly, the agreement between GSC and MC is the worst of all previously 
considered cases here. Even the values of I?o/c do notpmatch for small k because the core monomer is strongly affected 
by the very pronounced correlation hole effect in MCcJ. However, Vij between monomers within same arms and away 
from the core are naturally closer between GSC and MC, just as for the open chain in Fig. ||. Clearly, the divergence 
of the GSC and MC curves does not increase with k after fc '--^ 10 and the curves have rather similar overall shapes. 

It is interesting to analyse the structure of the collapsed globule now. Thus, in Fig. I the mean-squared distances 
7?ofc are plotted vs k for a flexible ring homopolymer in the globular state, [/'■''■' = 6. Overall shapes of the GSC and 
MC curves are quite similar, reflecting the compactness of the globule. The discrepancy between GSC and MC at 
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small k is present systematically. GSC generally overestimates the values of all X'ofci but the theory with a = 0.9 
manages to come to a nearly correct limit of the globule size for k* ^ 20 (see Ref. ^J), whereas the theory with a = 1.0 
is less accurate. The data for 37?,^ in the last block of Tab. || thus show a much better agreement than before, but 
the energy values in Tab. | have a significantly larger discrepancy between GSC and MC. The GSC theory noticeably 
underestimates the negative Lennard- Jones energy contribution between all pairs of monomers. This, however, is to 
be expected given that the shape of RDF from MCcJ has a very tall first liquid-like peak (see the inset of Fig. |^). 
Since the Lennard-Jones interaction is rather short-ranged, the first peak gives a predominant negative contribution 
to the mean energy in Eq. As the GSC theory has merely an effective smooth 'interpolating' Gaussoid in g^'^^ in 
that area (see the inset of Fig. , the resulting negative energy contribution is significantly smaller in such a theory. 

We can also remark that the scaling for the swelling exponent in the globule is correct in the GSC theory. Indeed, 
the maximal compression is reached when r/ ;$ 1, therefore by considering z/ij — > oo in Eqs. ( p3| , p^ , p6[ ) we indeed 
obtain, ~ (rf +rf' )2iV2/3. 

Next, we would like to look at the plots of the mean-squared radius of gyration, 37?.^, and of the mean energy, 
g{mt) a^CToss the coil-to-globule transition. These are depicted in Figs. |^ and|| respectively. The following points can 
be made. First of all, the shapes of these curves are quite similar for MC and GSC with both values of a. Secondly, 
the coil-to-globule transition is continuous in all three cases, with the energy slope changing noticeably at around the 
theta-point. Thirdly, the transition occurs at a somewhat higher value of the attraction constant J7'^°^ in the GSC 
theory than in the MC simulation. This can be explained by the underestimation of the Lennard-Jones attraction 
energy in the globule discussed above. Lastly, 3TZg of the globule is partly overestimated by the GSC method with 
a = 1, since the value C/'"' = 6 is much closer to the point of the coil-to-globule transition for the theory with a = 1 
than that with a = 0.9 or MC. 

Finally, to understand the iV-dependence of the coil-to-globule transition in Fig. ^ we present the plots of the 
specific energy slope, N-^dS'-"^^'^ / dU^^\ vs C/(°) for the flexible rings of different sizes. These curves nearly coincide 
in the repulsive coil region, starting to diverge from a value of U^^'^ ^ 1.2. The region of the transition, where 
the quantity N~^dE'^"^^^ / dU'^'^^ experiences the most dramatic drop, becomes increasingly narrower with increasing 
polymer size N . Moreover, since the magnitude of the overall change in the specific energy slope also increases with 
iV, the coil-to-globule transition becomes 'sbaijpfr' with N in the GSC theory, consistent with the MC simulation 
data and the transition being of second orderlljll3. Note also, that the theta-point, which we may define e.g. as the 
point of the maximal change in N~'^dE'^™'''^ / dU'^^\ shifts towards lower values of ?7*^°^ with increasing iV. 

V. CONCLUSION 

In this paper we have developed a version of the Gaussian self-consistent (GSC) technique which does not rely 
on the virial-type expansion of the Hamiltonian in terms of powers of the density, J dr p{r)^ . As a result, it is now 
possible to apply the new method to practically any polymer model involving conventional molecular interactions. 

Thus, we have been able to compare the mean spatial characteristics and the energy values between the results 
from the GSC theory and Monte Carlo (MC) simulation based on precisely the same model, which includes the 
harmonic bonded and the Lennard-Jones pair-wise interactions. This comparison has been performed for three types 
of macromolecular architectures of an isolated chain: a ring, an open polymer, and a star. We have done this also 
across the range of the coil-to-globule transition, as well as for different degrees of polymerisation and degrees of 
flexibility. Naturally, the GSC theory agrees with the second order nature of the coil-to-globule transition for flexible 
homopolymers. _ . 

Importantly, the speed of numeric computation is much fastera in the GSC method than in the equivalent MC 
simulation for obtaining the same data, particularly so for systems possessing an extra kinematic symmetry, such as 
for rings or stars. 

Overall, the agreement in the shapes of the curves and many of the particular numerical values of observables 
between GSC and MC is better than one could have anticipated given the relative simplicity of the GSC technique. 
Where any significant level of deviation does occur, it has been identified as either related to the correlation hole effect 
at smalLseparations, or to the stretched Gaussian behaviour at large separations, in the radial distribution function 
(RDF)EJ. 

In particular, for the repulsive coil, the energy and the mean-squared distances between near monomers along 
the chain are quite accurate in the GSC theory as compared to the MC data (with a typical deviation of several 
percent), but the distances between remote monomers, and hence the radius of gyration, are overestimated by the 
GSC method (with a typical deviation of a few dozen percent). This is a well known drawback of such a theory, 
related to the overestimation of the Flory swelling exponent for long chains, due to the fact that the RDF here does 
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not have a stretched exponential behaviour at large separations. On the contrary, for a rather stiff coil, the GSC 
theory underestimates the radius of gyration. This may be due to shortcomings of the hard sphere Carnahan Starling 
term. For the collapsed globule, on the other hand, the distances and the radius of gyration are quite accurate in the 
GSC theory, although the mean energy is less so because of the lack of a sharp liquid-like peak in the RDF. 

To make the agreement of GSC and MC better, one has to finally overcome the most restrictive feature of the 
method the Gaussian shape of RDF itself. One possible way of doing this is to take a linear superposition of 
the Gaussian trial functions, thereby permitting 'stretching' of the Gaussoid. This should be sufficient for curing the 
problem with the swelling exponent of the repulsive coil in the GSC theory. Work along these lines is currently in 
progress. The main difficulties in doing this are in the considerable mathematical complications when calculating the 
non-Gaussian conformational entropy of the chain, as well as in the added numerical complexity, since a radial mesh 
for RDF would have to be introduced. Nevertheless, we intend to resolve these issues and hope to present a more 
accurate, and at last a non-Gaussian self-consistent theory in the near future. 
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FIGURE CAPTIONS 



FIG. 1. Plots of the half-ring radial distribution functions 5o^j^/2('') ^ ^ units) from MC (solid thick lines) and the GSC 
theory with o = 0.9 (thin dashed lines) vs the radial separation r (in £ units) for the flexible ring homopolymers with A = 
and K = 1. The main paxt of the figure corresponds to the good solvent U^°^ = and the polymer size N = 300, whereas the 
inset corresponds to the poor solvent U^°^ = 6 and the polymer size TV = 200. 

FIG. 2. The mean-squared distances 2?ofc (in units) of an open flexible homopolymer with A'^ = 200, A = 0, and k = 1 
in the good solvent, U^°^ = 0, vs the chain index k. Here and below the solid thick lines correspond to the MC data, the solid 
thin lines to the GSC theory with a = 1.0, and the dashed lines to the GSC theory with a = 0.9. 

FIG. 3. The mean-squared distances Vok (in £^ units) of a fairly flexible, A = 1, ring homopolymer with N = 300 and 
K = 1 in the good solvent, U^°^ = 0, vs the chain index k. 

FIG. 4. The mean-squared distances I>ofc (in £^ units) of a semi-flexible, A = 5, ring homopolymer with N = 300 and k = 1 
in the good solvent, U^°^ = 0, vs the chain index k. 

FIG. 5. The mean squared distances Dok (in units) from the core monomer of a flexible, A = 0, homopolymer star with 
/ = 12 arms, k = 1, and the arm length N/f = 50 in the good solvent, (y'"' = 0, vs the chain index k. 

FIG. 6. The mean-squaxed distances 2?ofc (in £^ units) for the globule of a flexible, A = 0, ring homopolymer with N = 200, 
K = 1, and U^°^ = 6 vs the chain index k. 

FIG. 7. The mean squared radius of gyration 3TZ^ (in £^ units) of a flexible, A = 0, homopolymer ring with N = 150 and 
K = 1 vs the dimensionless degree of the Lennard-Jones attraction, 17'"', across the coil-to-globule transition. 

FIG. 8. The mean energy f (in ksT units) of a flexible, A = 0, homopolymer ring with N = 150 and k = 1 vs the 
degree of the Lennard-Jones attraction, f/'"-' , across the coil-to-globule transition. 

FIG. 9. Plots of the specific energy slope, N~^d£^''"*'' /dU^°\ (in ksT units) of flexible, A = 0, homopolymer rings with 
K = 1 vs the degree of the Lennard-Jones attraction, J7'°' , across the coil-to-globule transition for different polymer sizes 
N = 50, 100, 200, 300 (firom top to bottom). These are obtained from the GSC theory with a = 0.9. 
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TABLES 



TABLE L Comparison of the mean energy Vjaiues £;(*"*) (in ksT units) for different homopolymers from the MC simulation 
(second column) based on the data set of Ref.O and from the GSC theory with the parameter a — 0.9 (third column) and 
a = 1.0 (fourth column). The fifth and sixth columns contain the relative deviation 5 = {£'^'"^\GSC)/£{MC) - 1) x 100% for 
these two cases given as a percentage. The model parameters, which are suppressed within the tables, such as e.g. A, U^'^\ are 
equal to zero. 
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TABLE II. Comparison of the mean-squared radius of gyration values, 3TZ^ (in units) for different homopolymers 
from the MC simulation (second column) based on the data set of Ref.Ej and from the GSC theory with the parame- 
ter a = 0.9 (third column) and a — 1.0 (fourth column). The fifth and sixth columns contain the relative deviation 
5 = {TZl{GSC)/nl{MC) - 1) X 100% for these two cases given as a percentage. 
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